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O: 

rj . ABSTRACT 

We present an analysis of the relation between the masses of cluster- and group-sized 
i-C ' halos, extracted from ACDM cosmological N-body and hydrodynamic simulations, 

and their velocity dispersions, at different redshifts from z = 2 to z = 0. The main 
q aim of this analysis is to understand how the implementation of baryonic physics in 

?-h , simulations affects such relation, i.e. to what extent the use of the velocity dispersion 

C/3 ■ as a proxy for cluster mass determination is hampered by the imperfect knowledge of 

' the baryonic physics. In our analysis we use several sets of simulations with different 

physics implemented: one DM-only simulation, one simulation with non-radiative gas, 
^vq ' and two radiative simulations, one of which with feedback from Active Galactic Nuclei. 

| Velocity dispersions are determined using three different tracers, dark matter (DM 

■ hereafter) particles, subhalos, and galaxies. 

OO ' We confirm that DM particles trace a relation that is fully consistent with the 

\& , theoretical expectations based on the virial theorem, a v oc M a with a = 1/3, and with 

• previous results presented in the literature. On the other hand, subhalos and galaxies 

' trace steeper relations, with velocity dispersion scaling with mass with a > 1/3, and 

(~) with larger values of the normalization. Such relations imply that galaxies and subhalos 

(T) ■ have a ~ 10 per cent velocity bias relative to the DM particles, which can be cither 

t-H | positive or negative, depending on halo mass, redshift and physics implemented in the 

. simulation. 
i-H ' We explain these differences as due to dynamical processes, namely dynamical 

>^ \ friction and tidal disruption, acting on substructures and galaxies, but not on DM 

5— i ■ particles. These processes appear to be more or less effective, depending on the halo 

masses and the importance of baryon cooling, and may create a non-trivial dependence 
of the velocity bias and the C1D-M200 relation on the tracer, the halo mass and its 
redshift. 

These results are relevant in view of the application of velocity dispersion as a 
proxy for cluster masses in ongoing and future large redshift surveys. 

Key words: galaxies: clusters: general - galaxies: groups: general - galaxies: kine- 
matics and dynamics - galaxies: evolution - methods: numerical - cosmology: theory 



1 INTRODUCTION 

Galaxy clusters provide a powerful means of tracing the 
growth of cosmic structures and, ultimately, c onstraining 
cosmological parameters (e.g. lAllen et al.ll201ll ). A crucial 
aspect in the cosmological application of galaxy clusters con- 
cerns the reliability of mass estimates. Mass is not a di- 
rectly observable quantity but can be determined in several 



ways, e.g. by assuming the condition of e quilibrium of the 
intracluster plasma (e.g. lEttori et al]|2002T ) or galaxies (e.g. 
iKatgert et alj|2004 ) within the cluster potential well, or by 
measuring the gravitational lensing distortion of the images 
of background galaxies by the cluster gravitational field (e.g. 
lHoekstrdl2003h . 

These methods of mass measurement can only be ap- 
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plied to clusters for which high quality data are avail- 
able. When these are not available, it is still possible to 
infer cluster mass from other observed quantities, the so- 
called mass proxies, which are at the same time relatively 
easy to measure and ch aracterised by tight scaling relations 
with cluster mass (e.g. iKravtsov fc Borganil 12012 ). Exam- 
ples of such mass proxies are the total thermal content 
of th e intra-cluster plasma, measured from either X-ray 
(e.g., IKravtsov et al.|[2006l ; IStanek et ah| 2010l; iFabian et all 
20101) or Sunvaev-Ze l' dovich (e . g. lArnaud et all I20T0F 
Williamson et a.1.1 l201ll : iKav et all [2012! ) observations, the 
optical luminosi ty or richness traced by the cluster g alaxy 
population (e.g. iPopesso et al.ll2007l ; iRozo et al.ll2009|). and 
veloc i ty dispersion of member galaxies (e.g. Biviano et al.l 
l200d : ISaro et alj|2012j ). 

The use of velocity dispersion as a proxy for cluster 
mass is particularly interes ting in view of ongoing (BOSS, 
see, e.g. White et al.l l201ll ) and forthcoming (Euclid, see 
lLaureiis et al. 2011 ) large spectroscopic galaxy surveys. It 
is crucial to understand whether a cluster velocity disper- 
sion measured on its member galaxies is a reliable proxy for 
its mass. Calibration of such scaling relation can be based 
on detailed multi-wavelength observations of control sam- 
ples of galaxy clusters. On the other hand, detailed cosmo- 
logical simulations are quite useful to calibrate such scaling 
relations indepe ndently from possible observ ational system- 
atic effects (e.g. iBorgani fe Kravtsovl l201ll . and references 
therein) . 

The implementation of baryonic physics can play a fun- 
damental role in these analysis. In principle, since galaxies 
are nearly collisionless tracers of the gravitational potential, 
one expects velocity dispersion to be more robust than X- 
ray and SZ mass proxies against the effects induced by the 
presence of baryons and by their thermal history. 

Using a set of cluster-size d halos extra c ted f rom a 
ACDM cosmological simulation, iBiviano et al] i|2006h anal- 
ysed the reliability of the velocity dispersion as a mass proxy. 
They considered both DM particles and simulated galaxies 
as tracers of the gravitational potential of their host halo. 
They found that in typical observational situations, the use 
of the line-of-sight velocity dispersion, trios > allows a more 
precise cluster mass estimation than the use of the virial 
theorem. They used only one kind of simulation without ex- 
plor ing different baryon ic physics. 

lEvrard et al.l (|2008l ) analysed the 0200 — M200 El relation 

using several cosmological simulations, and showed that it 

1/3 

is close to the virial scaling relation 0200 oc M 2 q across 
a broad range of halo masses, redshifts, and cosmological 
models. When looking at simulated galaxies, some studies 
found that they show a significant, albeit small, velocity 
bias with respect t o DM particles jDiemand et all 



Faltenbacher et al. 2005 



Faltenbacher fc DiemandT 



2004 



2006 



Faltenbacher fc Mathews] 120071 ; lLau et all l20icf h These 



studies agreed that the amplitude of the velocity bias, i.e. the 
ratio between the velocity dispersions of simulated galaxies 



In the following we indicate with T200 the radius of the sphere 
drawn from the halo centre, which encloses a mean density of 
200 X p c (z), with p c (z) being the critical cosmic density at redshift 
z. The quantities with the subscript 200 have to be considered 
evaluated at, or within, r2oo- 



and DM particles, is not larger than « 10%, but disagreed 
on whether galaxies are positively (velocity bias > 1) or 
negatively (velocity bias < 1) biased. The disagreement is 
unlikely to come from resolution issues l|Evrard et al.ll2"008r ). 
Other effects are more important in affecting the value 
of the bias, such as t he distance from the cluster centre 
|Diemand et all |2004| ; iGill et al.l |2004). b arvon dissipation 
and redshift dependence ( Lau et al. 2010l ). 

The way simulated galaxies are selected also has an im- 
portant effect on the amount of velocity bias. By select- 
ing simulated galaxies in stellar mass, rather than in to- 
tal mass, the velocity bias is strongly reduced or even sup - 
pressed (|Faltenbacher fc Diemandl I2OO6I ; lLau et all 120101 ) . 



A similar effect is seen when the selection of simulated 
galaxies is based on their total mass at the moment 
of infall, which is found to be proportional to the stel- 
lar mass dFaltenbacher fc Diemandl l200d ; lLau et all l20ld : 
IWetzel fc White! l20ltf ). The proportional ty between to- 
tal and stellar mass of a galaxy is lost after the galaxy 
enters the cluster, because the DM halo is more easily 
strippe d than the stellar component by the cluste r tida l 
forces JDiemand et al.l |2004|; IBovlan-Kolchin et al.l 120081 ; 
lLauet all l201fj| ; IWetzel fc White! l2010h . Observational evi- 
dence for tidal strippi ng of cluster galaxies has been obtained 
from lensing studies (|Nataraian et al]|2002l ; [Limousin et al.l 
2007). Tidal stripping i s more effective for g alaxies mov- 
ing at lower velocities l|Diemand et al.l l2004h . When the 
mass removed from a simulated galaxy by tidal stripping is 
such that the galaxy mass drops below the resolution limit, 
the galaxy is effectively tidally disrupted. As the galaxies 
that are disrupted are preferentially those of smaller veloc- 
ities, the survivers will display on average a larger velocity 
dispersion than DM particles, i.e. a positive velocity bias 
|Faltenbacher fc Diemandll2006h . 

Another important process is dynami cal friction 
jChandrasekharl 1 19431 ; lEsquivel fc Fuchs! |2007| ). which re- 
moves energy from a galaxy orbit, bringing it closer to 
the cluster centre, and slowing down its velo city (e.g. 
IBovlan-Kolchin et allliooj ; IWetzel fc White!l20ld ). If a suf- 
ficient number of galaxies is slowed down by dynamical fric- 
tion and survive both tidal disruption and merging with the 
central galaxy, dynamical friction might cause a negative 
velocity bias in the cluster galaxy population. 

All the processes discussed so far alter the dynamics 
of tracers like galaxies, providing a source of uncertainty in 
the aforementioned relation linking mass and velocity dis- 
persion. With that comes the need of further investigations 
about this topic, also because of the different results found 
in the literature on the velocity bias of cluster galaxies. The 
aim of this paper is to characterise the relation between the 
velocity dispersion and the mass of simulated halos span- 
ning a wide mass range, from ~ 10 13 Mq to > 1O 15 M0, at 
different redshifts (from z = 2 to z = 0), and using differ- 
ent tracers, DM particles, subhalos, and galaxies, in order 
to understand how reliable is the velocity dispersion as a 
proxy for cluster masses. Simulations with different physics 
implemented are used, in order to understand how differ- 
ent physical processes affect the structures and hence the 
dynamics of tracers. 

In this paper we do not consider observational bi- 
ases, such as projec t ion effects and p resence of interlopers 
l|Biviano et al.ll2006l ; ISaro et al.ll2012T l. but we focus on the 
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effects due to the physics and the implementation of bary- 
onic physics in the simulations. We find that such implemen- 
tation affects the kinematics of the systems. The analysis of 
observational effects must therefore be based on simulations 
where baryonic physics is taken into account. 

This paper is structured as follows. In Sect. [2] we de- 
scribe the simulations used for this work and define the 
samples used in our analyses. In Sect. we determine the 
relation between mass and velocity dispersion, and how it 
depends on redshift and on the different types of simula- 
tions. In Sect. 13.11 we quantify the scatter and its nature 
(statistical or intrinsic). In Sect. 13.21 we describe the velocity 
bias of subhalos and galaxies with respect to the underlying 
diffuse component of DM particles. In Sect. [331 we look for 
a signature of the different dynamical processes which are 
at work in galaxy systems, on the velocity distributions of 
the different tracers of the gravitational potential. Finally in 
Sect. [4] we discuss our results and present our conclusions. 



2 SIMULATIONS 
2.1 Initial conditions 

Our samples of cluster-sized and group-sized halos are 
obtained from 29 Lagrangian regions, centred around as 
many massive halos identified within a large-volume, low- 
resolution N-body cosmological simulation, resimulate d with 
higher resolution. We refer to iBonafede et al.l |201lft for a 
more detailed description of the set of initial conditions used 
to generate samples of simulated clusters used for our anal- 
ysis. 

The parent Dark Matter (DM) simulation followed 
1024 3 DM particles within a box having a comoving side 
of 1 h -1 Gpc, with h the Hubble constant in units of 
100 km s" 1 Mpc -1 . The cosmological model assumed is 
a flat ACDM one, with Q m = 0.24 for the matter den- 
sity parameter, fibar = 0.04 for the contribution of baryons, 
Ho = 72 km s _1 Mpc -1 for the present-day Hubble constant, 
n s = 0.96 for the primordial spectral index and as = 0.8 
for the normalisation of the power spectrum. Within each 
Lagrangian region we increased the mass resolution and 
added the relevant high-frequency modes of the power spec- 
trum, following the zo omed initial condition (ZIC) technique 
|Tormen et all 1 1997ft . Outside these regions, particles of 
mass increasing with distance from the target halo are used, 
so that the computational effort is concentrated on the re- 
gion of interest, while a correct description of the large scale 
tidal field is preserved. Each high-resolution Lagrangian re- 
gion is shaped in such a way that no low-resolution particle 
contaminates the central zoomed-in halo at z = at least 
out to 5 virial radii. As a result, each region is sufficiently 
large to contain more than one interesting halo with no con- 
taminants within its virial radius. 

Initial conditions have been first generated both for 
DM-only simulations. The mass of DM particles in the 
zoomed-in regions is hidm = 10 9 h _1 M Q . Henceforth we 
refer to these simulation as DM-only. Initial conditions for 
hydrodynamical simulations have been generated only in the 
low-resolution version, by splitting each particle within the 
high-resolution region into two, one representing DM and 
another representing the gas component, with a mass ratio 



such to reproduce the cosmic baryon fraction. The mass of 
each DM particle is then hidm = 8.47 • 10 s h _1 M and the 
initial mass of each gas particle is m gas = 1.53 • 10 8 h _1 M0. 

2.2 The simulation models 

All the simulations have been carried out with the TreeePM- 
SPH GADGET-3 code, a more efficie nt version of the previ- 
ous GADGET-2 code |SpringelH2005ft . As for the computa- 
tion of the gravitational force, the Plummer-equivalent soft- 
ening length is fixed to e = 5 h _1 kpc in physical units below 
z = 2, while being kept fixed in comoving units at higher 
redshift. 

Besides the DM-only simulation, we also carried out a 
set of non-radiative hydrodynamic simulations (NR here- 
after) and two sets of radiative simulations, based on differ- 
ent models for the release of energy feedback. 

A first set of radiative simulations includes star for- 
mation and the effect of feedback triggered by supernova 
(SN) explosions (CSF set hereafter). Radiative cooling rates 
are computed by follo wing the same procedure presented by 
IWiersma et al.l |2009r i. We account for the presence of the 
cosmic microwave background (CMB) and for the model of 
UV/X-ray backgro und radiation from qua sars and galax- 
ies, as computed bv lHaardt fc Madaul (|200lft . The contribu- 
tions to cooling from each one of eleven elements (H, He, 
C, N, O, Ne, Mg, Si, S, Ca, Fe) have been pre-computed 
using the publicly av ailable CLOUDY photo-ionisation code 
ijFerland et al.l [1998) for an optically thin gas in (photo- 
ionisation) equilibrium. Gas particles above a given thresh- 
old density are treated as multiphase, so as to provide a sub- 
resolution description of the interstell ar medium, according 
to the model originally described by ISpringel fe Hernquistl 
(|2003ft . We also include a description of metal produc- 
tion from chemical enrichment contribute d by SN-II, SN- 
Ia an d low and intermediate mass stars jTornatore et al.l 
2007). Stars of different mass, distributed according to a 
Chabrier IMF (|Chabrierl 12003ft . release metals over the 
time-scale determined by the corresponding mass-d ependent 
life-times (taken from iPadovani fc Matteuccilll993h . Kinetic 
feedback contributed by S N-II is implemented ac c ording to 
the scheme introduced by ISpringel fc Hernquistl (|2003ft : a 
multi-phase star particle is assigned a probability to be 
uploaded in galactic outflows, which is proportional to its 
star formation rate. In the CSF simulation set we assume 
v w = 500 km s -1 for the wind velocity. 

Another set of radiative simulations is carried out by in- 
cluding the same physical processes as in the CSF case, with 
a lower wind velocity of v w = 350kms _1 , but also includ- 
ing the effect of AGN feedback (AGN set, hereafter). In the 
model for AGN feedback, released energy results from gas 
accretion onto supermassive black holes (BH). This model 
introduces so me modifications with respect to that originally 
presented by ISpringel et all (|2005ft (SMH) and will be de- 
scribed in detail by Dolag et al. (2012, in preparation). BHs 
are described as sink particles, which grow their mass by gas 
accretion and merging with other BHs. Gas accretion pro- 
ceeds at a Bondi rate, while being Eddington-limited. Radi- 
ated energy corresponds to a fraction of the rest-mass energy 
of the accreted gas. This fraction is determined by the radi- 
ation efficiency parameter e r = 0.1. The BH mass is corre- 
spondingly decreased by this amount. A fraction of this radi- 
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ated energy is thermally coupled to the surrounding gas. We 
use e/ = 0.1 for this feedback efficiency, which increases to 
e/ = 0.4 whenever accretion enters in the quiescent "radio" 
mode and takes place at a rate smaller than one-hundredth 
of th e Eddington limit (e.g. lSiiacki et af1l2007l ; [Fabian et alj 
l2010h . 



2.3 The samples of simulated clusters 

The identification of clusters proceeds by using a catalogue 
of FoF groups as a starting point. The S UBFIND algorithm 
(|Springel et al J l200ll ; iDolag et af]|2009l ) is used to identify 
the main halo, whose centre corresponds to the position of 
the most bound DM particle, and substructures within each 
FoF group. In the following, we will name "galaxies" the 
bound stellar structures hosted by the subhalos identified by 
SUBFIND in the radiative CSF and AGN hydrodynamical 
simulations. 

In this work we consider all the main halos with M200 > 
10 13 Mq from z = to z — 2, which contain no low- 
resolution particles within r2oo- Among these cluster-sized 
and group-sized halos, we only retain those with at least 
five subhalos more massive than 10 M© within r2oo- The 
number of selected halos varies at different redshifts and in 
different simulation sets, from a minimum of 54 to a maxi- 
mum of 308. 

The subhalos that we consider in our analysis are se- 
lected to be more massive than 10 M©, which corresponds 
to 72 particles in the DM simulations. The galaxies we con- 
sider in our analysis (in the CSF and AGN sets) are selected 
to have a stellar mass ^ 3 x 10 9 M©. By choosing this lower 
limit in stellar mass, we retain all subhalos more massive 
than 10 1 M© and include many others with smaller masses. 
As a consequence, there are more halos with 5 galaxies 
than with ^ 5 subhalos (within r2oo)- Note that the effects 
of the AGN feedback is negligible for galaxies with stellar 
masses below the chosen limit since they are generally hosted 
within halos where BH particles have never been seeded. 



3 THE VELOCITY DISPERSION 
RELATION 



MASS 



Given a tracer of the gravitational potential of a halo (DM 
particles, subhalos, galaxies) it is possible to write a re- 
lation between halo mass and velocity dispersion of the 
tracer, based on (i) the definition of circular velocity at 
r200, V200 = 10 [G h(z) M200] , and (ii) the relation be- 
tween (T200 and V2oo- A relation between velocity dispersion 
and mass can be derived analytically once the form of the 
mass density pr ofile and of the velocity an isotropy profile are 
given (see e.g. Mauduit fc Mamonll2007l . and the erratum 
iMauduit fc MamorJ|2009l). 

Following lEvrard etal] | |200St ). we use the one- 
dimensional velocity dispersion trip = 0200/ V3. Using the 
relations provided bv lLokas fc Mamonl (|200ll . eqs. 22-26) we 
calculate the ratio between <tid and V200 for NFW mass pro- 
files of concentrations c = 3 and 10, and for (constant) veloc- 
ity anisotropiefl /3 = 1 — (at/a r ) 2 — and 0.5. These values 



of c and /3 are typical of group- and c l uster-sized halos (e.g . 
iGao et alj I200&) ; IWoitak et all l20Qg| ; iMamon et all |201ofi . 

We find ctid/v2oo = 0.64,0.69 for /3 = and <tid/v2<x> = 
0.68, 0.70 for j3 = 0.5; the larger values are for c = 10. The 
average value we calculate using DM particles for all the ha- 
los selected in our analysis, <tid/v2oo = 0.68, lies within the 
same range. 

Using the range ctid/v2O0 = 0.64-0.70 just found, and 
the definition of V200 , we find 



ClD 

km s~ 



= Ai 



h(z) M200 



10 15 M 







(1) 



with Am = 1040-1140 and a = 1/3. Given that the real, 
simulated halos are not perfect NFW spheres, and that a 
halo concentration depends on its mass and redshift (e.g. 
iDuffv et al.ll2008l : IGao et al J [20081 ). the real values of Aid 
and a can be different, and need to be evaluated from the 
data. Moreover, DM particles, subhalos, and galaxies do not 
necessarily obey the same (T1D-M200 relation. We therefore 
evaluate for each cluster of the samples described in Sect. 
I2.3l the values of uid of DM particles, subhalos, and (for the 
CSF and AGN simulations) galaxies. The use of simulations 
with different baryonic physics implemented allows us to un- 
derstand how baryons and different feedback models modify 
the scaling relation between velocity dispersion and mass. 
The (7 id values o f DM particles are obtained using the bi- 
weight estimator (|Beers et a.l.lll990h . when at l east 15 data 
pints are available. Otherwise, as suggested bv lBeers et al.l 
1990). we use the classical standard deviation. The confi- 
de nce intervals for th e o\ld values are obtained using eq. (16) 
in lBeers et all ljl990h . 

We then perform a linear fit to the log(uiD) vs. 
log[h(z) M200] values, with am in units of km s _1 and M200 
in units of 10 h _1 M©, for each simulation set, at several 
redshifts. The fits were performed with the IDL procedure 
linf it, inversely weighting the data by the uncertainties in 
the values of <tid. The results of these fits are the values 
of the parameters Aid and a of eq. ^Q. In Fig. [T] and [2] 
we show examples of these fits for the AGN simulations at 
redshift and 1.26, respectively. 

In Fig. [3] we show the best fitting values of Aid and a 
as a function of redshift, for different tracers in the AGN 
simulation, as an example. The slope a is confirmed to be 
consistent with the theoretically expected value a = 1/3. On 
the other hand, its value is significantly larger when using 
subhalos or galaxies as tracers. In any case the values of a 
and Am do not generally show a significant dependence on 
z. Only in two cases (subhalos in the DM simulation and 
galaxies in the AGN simulation - see Appendix [A} we do 
find (marginally) significant correlations between a and z, 
mostly driven by the point at z = 2. Even in these cases, 
a changes very little with redshift, ~ 4% for the galaxies in 
the AGN simulation from z=1.5toz = 0. In fact, a model 
where a is constant with z provides an acceptable fit (in a 
X 2 sense) to all cases (also those not shown in Fig.[3J). Since 
the variation of Aid and a with z is not significant, we take 
the weighted averages of their values over all redshifts to 
characterise the ctid — M200 relations of the different types 
of simulations and tracers (see Table [1] and Fig. [4}. 



The quantities at = [(o^ + cr^)/2] 1 / 2 and <r r are, respectively, 



the ID tangential and radial component of the 3D velocity dis- 
persion. 
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Figure 1. Top panels: velocity dispersion a^D (km s 1 ) as a function of halo mass h(z) M 20 o (10 15 M ), for DM particles (left panels), 
subhalos (central panels), and galaxies (right panels), at z=0, in the AGN simulation sets. The dashed green line represents the theoret- 
ically expected relation tr^ = 1090 ■ (h(z) M2uo) 1/ ' 3 ; the solid line in each panel represents the best-fit relation. Bottom panels: y-axis 
residuals of the DM particles (left), subhalos (centre), and galaxies (right), from the DM best-fit relation. 
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Figure 3. Average values of a and Aid (see eq. [TJ for halos in 
AGN simulations as a function of redshift, for different tracers, 
DM particle (black dots), subhalos (red diamonds), and galax- 
ies (blue triangles). The dashed green line in the top panel is 
the theoretically expected value a = 1/3, the horizontal shaded 
area in the bottom panel indicates the expected theoretical range 
Aid = 1040-1140 (see text). The data of galaxies have been 
shifted by —0.05 in z for the sake of clarity. 





A 1D (km s" 1 ) 


a 


DM 


1094 ±3.7 


0.334 ± 0.0014 


NR. 


1102 ±3.4 


0.336 ± 0.0021 


CSF 


1081 ±4.1 


0.329 ± 0.0021 


AGN 


1095 ± 4.4 


0.336 ± 0.0015 


DM sub 


1244 ± 4.7 


0.361 ±0.0027 


NR sub 


1259 ±5.3 


0.364 ± 0.0030 


CSF sub 


1166 ±5.1 


0.362 ±0.0018 


AGN sub 


1199 ±5.2 


0.365 ±0.0017 


CSF gal 


1165 ±6.7 


0.355 ±0.0025 


AGN gal 


1177 ±4.2 


0.364 ± 0.0021 



Figure 4. Values of a and Aid (see eq. [TJ averaged over all 
redshifts for different tracers - DM particle (black dots), subhalos 
(red diamonds), and galaxies (blue triangles) - for the different 
types of simulations. The dashed green line in the top panel is 
the theoretically expected value a = 1/3, the horizontal shaded 
area in the bottom panel indicates the expected theoretical range 
Aid = 1040-1140 (see the text). 



(|200ST ) and co ming from the analys is of different simulations 
(see Table 4 in lEvrard et all 120081), most of which are DM- 
only. The comparison with Lau et all (|201p| ) is less straight- 
forward, since they did find a z-dependence of the a pa- 
rameter for DM particles. Taking an error-weighted aver- 
age of the values they found at different redshifts, we ob- 
tain A 1D = 1103 ± 2 km s~\ a = 0.325 ± 0.011 for their 
non-radiative simulations, and Aid = 1160 ± 9 km s~\ 
a = .296 ± 0.012 for their radiative simulation. lLau et all 
(|20ld )'s 1T1D-M200 relations are therefore flatter than ours, 
particularly so for the radiative simulations. We discuss 
these differences in Sect. [4] 



Table 1. Weighted average values over all redshifts of the param- 
eters Am an d a defining the <J\d — M200 relation (see eq. JTJ) for 
different simulation sets, for DM particles (top table), subhalos 
(middle table), and galaxies (bottom table). 

In Fig. [4] we show the dependence of the parameters a 
and A\o on the physical processes included in the simula- 
tions. When considering DM particles as tracers, the Aid 
values are well within the theoretically expected range, and 
the a values are close to the virial expectation 1/3, regard- 
less of the baryonic physics implemented in the simulations. 
When considering subhalos or galaxies as tracers, the ctid- 
M200 relations are significantly steeper (a > 1/3) than for 
DM particles and than expected theoretically. Furthermore, 
while the slope is nearly the same for all the simulation sets, 
with a ~ 0.36, the same does not hold for the normaliza- 
tion Aid, this value being smaller when cooling and star 
formation are included. In each simulation set, the a and 
Aid values are highest for the subhalos, and lowest for the 
DM particles, with the values for the galaxies being closer 
to those of the subhalos. 

Both the Am and the a values we find fo r the DM par- 
ticles are in the range of the values listed bv lEvrard et al.l 



3.1 Scatter 

An analysis of the scatter of the scaling relations is quite 
important fo r cosmological applicati ons. It has been pointed 
out (see, e.g. iMortonson et al.ll201ll ) that the so called Ed- 
dington bias causes the mass function to shift toward higher 
values of mass if the scatter in the scaling relation between 
mass and mass proxy is not properly taken into account. Fur- 
thermore we want to understand the nature of such scatter, 
that is whether it is intrinsic or mainly due to statistical 
uncertainties. 

The best fit relation eq. {JJ has been subtracted from 
the values of the velocity dispersion of the clusters for each 
simulation at each redshift. The result for the AGN set at 
z=0 is shown in Fig. [5] for the three trac ers. The errors ar e 
associated to the biweight average value jBeers et al.lll990l ). 
In this figure the distribution of the residuals is also shown, 
along with the best fit gaussian curves and the reduced \ 2 
values of the fits. The residuals appear to be normally dis- 
tributed, substructures and galaxies having a wider distri- 
bution. 

The scatter appears to be independent of cluster mass, 
as well as of redshift and the type of simulation as shown in 
top panels of Fig. [6] the only difference being in the tracer 
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Figure 6. Left panel: Fractional scatter as a function of redshift 
in the AGN set. Black solid line refers to the value computed 
by using DM particles, red dashed and blue dash-dotted by using 
substructures and galaxies, respectively. The thinner lines are the 
confidence intervals. Right panel: Fractional scatter as a function 
of the type of simulation at z=0. Lines and colours are the same 
as in the left panel. 

used to build the (T1D-M200 relation. When using DM par- 
ticles the scatter is around 5%, while using substructures or 
galaxies it is around 12%. 

In order to understand the nature of the scatter, that 
is whether it is intrinsic or statistica l, we have tried t o com - 
pute the intrinsic scatter following IWilliams et~aH (|2010t ). 
Performing a linear fit of the logcrirj — logM2oo relation, 
the intrinsic scatter is considered as a parameter that mini- 
mizes x 2 = E,1l" [Vi ~(a-Xi + 6)] 2 /[ £ L + ° 2 ' e *,i + ^nt], 
where y — a ■ x + b is the linear relation, e x ,i and e Vf i are 
the uncertainties on the x and y quantities and ai„t is the 
intrinsic scatter. In order to estimate the value of the intrin- 
sic scatter and its uncertainty, we performed 1000 bootstrap 
resamplings of the couples of values (M200, 0"id), each time 
computing the intrinsic scatter estimate. In a first time, we 
have evaluated the intrinsic scatter in 4 bins of mass, but 
we found no mass dependence. Therefore we have evalu- 
ated it using all the data. The values of velocity dispersion 
evaluated using DM particles have been obtained using a 
huge number of objects, hence the statistical uncertainty is 
very small. The 5% scatter is therefore entirely intrinsic. On 
the other hand velocity dispersions estimates using substruc- 
tures and galaxies are based on relatively small numbers of 
objects, and statistical uncertainties dominate the scatter. 
The resulting intrinsic scatter for these tracers turns out to 
be quite small and consistent with the value found for DM 
particles. 

3.2 Velocity bias 

The difference between the CT1D-M200 relations established 
for DM particles, on one side, and subhalos and galaxies, 
on the other side implies that subhalos and galaxies have a 
higher velocity dispersion than DM particles in high mass 
halos. Since the relation is steeper for subhalos and galaxies 
than for DM particles, the opposite may occur in halos of 
low masses. To examine how am changes in halos of diflfer- 



Figure 7. Average values of the velocity bias for subhalos, as a 
function of redshift, for halos in the low-mass sample (left panel) 
and in the high- mass sample (right panel). Dotted line refers to 
substructures in the DM set, short black dashed line refers to 
substructures in the NR set, the other lines refer to substructures 
of the CSF and AGN sets, as indicated in the legend. Error-bars 
on the bias values are not shown for the sake of clarity. A typical 
error-bar is shown at the bottom-right of the right panel. The 
other points with error bars are the z = values from lLau et all 
(2010). The legend in the right panel identifies the different sym- 
bols as representative o f the non-radiativ e simulations (NR), and 
the CSF simulations of lLau et all l l2010h . 
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Figure 8. Same as Fig. [7] but for galaxies rather than substruc- 
tures. 



ent masses when using different implementations of baryonic 
physics, here we analyse the velocity bias, i.e. the ratio be- 
tween the velocity dispersions of subhalos (galaxies) and DM 
particles, b v = a sub /a D M (b v = a-gai/o"DM, respectively), 
for two subsamples of halos, one with masses h(z)M2oo < 
1O 14 M , and the other with masses h(z)M 20 o > 3 • 10 14 M o 
(the low- and high-mass samples hereafter). The average b v 
values as a function of redshift are shown in Fig. [7] and |SJ 
for subhalos and galaxies in different simulation sets, sepa- 
rately for the halos in the low-mass sample (left panels) and 
for the halos in the high-mass sample (right panels). 

The bias is on average higher for the high-mass sample 
halos, and at lower redshifts. A remarkable change in the 
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bias appears when introducing cooling and star formation. 
In fact the biases are higher for substructures in the DM 
and in non-radiative simulations (NR) than in the radia- 
tive ones (CSF, AGN), both when using substructures and 
galaxies, reflecting again the importance of the cooling and 
the feedback processes in the dynamics. 

In Fig. [7] a nd [5] we also show the values found by 
lLau et all (|201Ch . At z = 0, the sample of lLau et all (|2010T ) 
and ours have similar mass distributions and there is a good 
agreement (within the error bars) between their b v values 
and ours, separately for the different types of simulations 
and for the differe nt tracers. The c omparison between our 
data and those in lLau et aL (2010) at z > is not shown 
because it is not straightforward: Lau et aTl (|2010l ) follow 
the same halos at different redshifts by always considering 
the most massive progenitors of the halos selected at z = 0, 
whereas we consider all halos above a given mass cut at any 
given redshift. Taken at face value, the biases they find at 
higher z are smaller than those we find, and this difference 
might be related to a strong overcooling present in their 
simulations, making their subhalos very resistant to tidal 
disruption. 



3.3 Dynamical processes in halos 

The above results show how the relative kinetic energy con- 
tent in subhalos (or galaxies) and the diffuse DM compo- 
nent varies with redshift, halo mass, and the type of simula- 
tions. Here we provide an interpretation of these results in 
terms of the dynamical processes that are effective in clusters 
and groups, i.e. dynamical friction, merging with the central 



galaxy, tidal disruption, and accretion from the surrounding 
large-scale structure, and how these processes depend on the 
physics implemented in the simulations. 

A better understanding of the physical cause of the dif- 
ferences in the (T1D-M200 relation and of the velocity biases 
analysed in the two previous sections can be obtained by ex- 
amining the probability distribution function for the modu- 
lus of the velocity of the different tracers. These distributions 
are shown in Fig. [9] for the CSF and AGN simulations, at z 
= and 1.26, separately for the low- and high-mass samples. 

The velocity distribution of DM particles (solid black 
lines in Fig. [9} appears to be single-peaked, but only for the 
low-mass sample. The DM particle velocity distribution for 
the high-mass sample is flat-topped and sometimes double- 
peaked. These features, present in all the different types of 
simulations, appear to be more pronounced for the veloc- 
ity distributions of subhalos (dashed red lines) and galaxies 
(dotted blue lines) , and depend on the physics implemented 
in the simulations. The difference in the velocity distribu- 
tions of DM particles, subhalos, and galaxies is at the origin 
of the velocity biases and the differences in their CT1D-M200 
relations. Hereafter we interpret the effect of the dynamical 
processes that shape these velocity distributions. 

A flat-topped or double-peaked velocity distribution is 
the signature of ongoing mergers and infall of matter into 
the halos, which tends to populate t he high-velocity part 
of the velocity di s tributions (see, e.g. | Diemand et alj |2004| ; 
IWang et all 120051 ; IWetzell l201ll ; iRibeiro et alj|201ll ). With 
time, the velocity distribution relaxes via phase-mixing and 
dynamical relaxation and approaches a single Maxwellian 
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(e.g. iFaltenbacher fc Diemandll2006l ; lLapi fc Cavalierell201ll . 

and references therein). 

The relative importance of the relaxed and unrelaxed 
velocity distributions depends on how strong the matter 
infall rate is. Higher-mass halos are dynamically younger, 
and undergo significant matter i nfall until more recent 
times than lower-mass halo s (e.g. lLapi fc Cavalierel 120091 : 
IFaltenbacher fc White! l2010h . Hence in higher-mass halos 
the unrelaxed, high-velocity component is more important 
than in lower-mass clusters, as we see in Fig. [9] 

The velocity distributions of subhalo and galaxies is dif- 
ferent from those of DM particles, because of additional 
dynamical processes that leave the DM particle distribu- 
tions unaffected. One is the dynamical friction, leading to 
a decrease in the orbital energy, hence to an approach to 
the cluster centre and a decrease in velocity, eventually 
followed by merger with t h e central cluster galax y (e.g. 
iBovlan-Kolchin et alj 120081 : IWetzel fc White! |2010T ). The 
other process is tidal disruption, caused by the integral of 
tidal interactions leading to mass losses (e.g. iDiemand et al.l 
I2004T ). These two processes are related. On one hand, dy- 
namical friction becomes ineffective when tidal mass losses 



become important (e.g. IFaltenbacher fc Mathewdl20"07T ). be- 
cause dynamical friction is more effective for more massive 
objects. On the other hand, tidal mass losses are enhanced 
by dynamical friction, since they are more effective in slow- 
moving s ubhalos (and galaxies ) than their fast-moving coun- 
terparts (|Diemand et alj|2004l ). Hence dynamical friction is 
likely to be more effective at the first or bit of a subhalo (or 
galaxy) |Faltenbacher fc "M athews 2007), while tidal disrup- 
tion may take several orbits. 

Dynamical friction tends to increase the low-velocity 
tail of the distribution. A possible signature of this effect is 
visible in the velocity distributions of subhalos and galaxies 
at high-z in Fig. [5] On the other hand, the removal of the 
slowest subhalos (or galaxies) by tidal disruption or merger 
with the central galaxy tends to decrease the low-velocity 
tail of the distribution. This is likely to occur in the re- 
laxed component of the global velocity distribution, since 
subhalos (and galaxies) in the relaxed component have spent 
more time orbiting their parent halos than the recently in- 
fallen, unrelaxed population. The resulting velocity distribu- 
tion will then appear to be double peaked, one low-velocity 
peak being due to the relaxed component, depopulated by 
tidal stripping, and another high-velocity peak due to the in- 
falling population, that has not yet experienced significant 
tidal mass losses. This is particularly evident in higher-mass 
halos, in which the infall rate is higher (at given z) than 
in lower-mass halos. The resulting asymmetric velocity dis- 
tribution is visible in Fig. [5] (dashed red and dotted blue 
lines). 

These processes are dependent on the physical processes 
and feedback implemented in the simulations. In fact in Fig. 
[9]one can note that the velocity distributions of galaxies and 
subhalos in the CSF simulation sets are both closer to the 
velocity distributions of DM particles, than the correspond- 
ing ones in the AGN simulation sets. This occurs because 
baryon cooling tends to make galaxi es and subhalos more 
resistant against tidal disruption (e.g. IWeinberg et al. I l2008l ; 
lLau et ail 201(j ). However cooling efficiency is reduced in 
simulations including AGN feedback thus making halos less 
compact and galaxies more fragile than in the over-cooled 
CSF simulations. 

Yet another difference is visible in Fig. [9] namely the 
velocity distribution for galaxies in the AGN sets is slightly 
closer to the distribution of DM particles than that of sub- 
halos. This is due to the different mass limit we use to select 
galaxies and subhalos (see Sect. [21 • In our samples, subha- 
los selected by means of total bound mass are on average 
objects of higher mass than subhalos selected by means of 
stellar mass, hence they are subject to stronger dynamical 
friction, slowing them down and making them more easily 
su bject to tidal disruption. Our res ult is reminiscent of that 
of IFaltenbacher fc Diemandl (|2006f ) who found that the ve- 
locity distribution of galaxies is quite similar to that of DM 
particles when selecting objects by their mass at accretion, 
i.e. before tidal stripping operates, as pointed out also by 
lLau et all (|2010h . 

The processes so far discussed provide a frame for a 
better understanding of Fig. [7] and [S] At fixed galaxy mass, 
dynamical friction is more efficient in low-mass clusters and 
at high redshifts, before tidal stripping decreases galaxy 
masses. Therefore dynamical friction tends to create a ve- 
locity bias < 1 in low mass, high-z systems, while leaving al- 
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Figure 10. Velocity anisotropy as a function of the cluster cen- 
ter, for the AGN set. Top panels refer to z=0, bottom panels to 
z=1.26. Panels on the left are relative to the low-mass sample, 
panels on the right to the high-mass sample. Red points refer to 
the results obtained using subhalos, while blue diamonds using 
galaxies. The four bins are built in such a way that within each 
bin there is the same number of objects. Galaxies diamonds have 
been shifted by -0.05 in rcdshift for the sake of clarity. The black 
solid line refers to the results obtained using DM particles. 



most unaffected high mass clusters. Tidal stripping is more 
efficient on slow moving galaxies, which are stripped and 
eventually disrupted (or removed from a mass-limited sam- 
ple), and operates at all times. This causes an increase of 
the bias with time. This has the effect of erasing the initial 
dynamical-friction bias in low-mass clusters and creating a 
bias > 1 in high-mass clusters as we approach z = 0. 

An important role in these processes may also be played 
by galaxy orbits. In Fig.QIJ]we show the velocity anisotropy 
profiles for the different tracers, in low- and high-mass sys- 
tems, for z=0 and 1.26. For lack of sufficient number of sub- 
structures and galaxies in each system, the profiles for these 
tracers are computed for stacks of all systems, where sub- 
structure and galaxy velocities have been scaled by each sys- 
tem V200 before stacking. With this procedure, richer clus- 
ters weigh more in the final profile. For the DM particles we 
are not limited by poor statistics, so we derive the anisotropy 
profiles individually for each halo, and then take an average. 
For consistency with what was done for the substructures 
and galaxies, the average is weighted by the number of sub- 
structures present in each cluster. Fig. [10] shows that orbits 
are more radially anisotropic i n high-mass than in low-mass 
systems (as already shown by IWetzelll201 il l . Understanding 
the reason for this difference (and for the redshift evolution 
clearly visible in the same figure) is beyond the scope of this 
paper, although we suspect that it might be related to the 
younger dynamical age of higher m ass clusters (as suggested 
by, e.g.. iBiviano fc Poggiantll2009l . and references therein). 
What is relevant in this context is that galaxies on more 
radial orbits reach closer to the cluster centre and therefore 
suffer from stronger tidal stripping effects. The different or- 
bital anisotropy of galaxies contribute to create a higher 
velocity bias in high-mass relative to low-mass clusters. 



4 DISCUSSION AND CONCLUSIONS 

We determined the (T1D-M200 relation for DM particles, sub- 
halos, and galaxies in cluster- and group-sized halos ex- 
tracted from ACDM cosmological simulations. We analysed 
four sets of simulations carried out for the same halos: one 
DM-only, one with non-radiative gas, and another two with 
gas cooling, star formation and galactic ejecta triggered by 
SN winds, one of the two also including the effect of AGN 
feedback. 

The main results of our analysis can be summarised as 
follows. 

• We confirm that the o\d — M200 relations for the DM 
particles are consistent with the theoretical expectation from 
the virial relation, assuming NFW halo mass profiles, with 
reasonable values of concentration and velocity anisotropy. 

• The intercepts at 1O 15 M0 and slopes of the logarithmic 
(T1D-M200 relations derived using subhalos and galaxies as 
tracers of the potential are always higher than those derived 
using DM particles. We do not find a significant dependence 
of the CT1D-M200 relations on redshift, but we do find a depen- 
dence on the kind of simulation, the radiative ones having a 
higher value of the normalization. 

• The CT1D-M200 relati ons for the DM parti cles are consis- 
tent with those found bv lEvrard et all (|200ST ). On the other 
hand the relations we find for all the tracers are steeper than 
those derived by lLau et ail (|201(t ). This difference might 
be caused by th e narrower halo mass range explored by 
lLau et al.l (|2010h . 

• The differences in the (T1D-M200 relations for the differ- 
ent tracers of the halo gravitational potential and for the dif- 
ferent physics implemented are due to dynamical processes 
taking place in the halos. In fact dynamical friction and tidal 
disruption act on substructures and galaxies but not on DM 
particles. Dynamical friction slows down a substructure or 
a galaxy before it suffers mass loss due to tidal stripping. 
Dynamical friction is therefore more efficient at high-z. It 
is also more efficient in lower-mass clusters for given galaxy 
mass. Tidal stripping, on the other hand, is more efficient 
in higher-mass clusters, where galaxies move on more ra- 
dial orbits (hence with smaller pericenter radii). As a result, 
velocity biases are created in the substructure and galaxy 
populations, relative to the DM particles, and these biases 
are ^ 1 for low-mass systems and ^ 1 for high-mass systems, 
and increasing with time, leading to the observed differences 
in the CT1D-M200 relations of different tracers. 

• In order to correctly reproduce such processes, a de- 
tailed implementation of the baryonic physics must be used 
in the simulations. In fact the presence of baryonic mat- 
ter makes halos more resistant against tidal disruption (e.g. 
IWeinberg et al.ll2008f h In this way in simulations with cool- 
ing and star formation there is a higher fraction of survivors 
among the slow moving subhalos, reflecting in a lower value 
of the normalization on the CT1D-M200 relation. 

• We analysed the scatter in the am — M200 relation, 
finding a value of around 5% for DM particles and 12% for 
substructures and galaxies. The intrinsic scatter, after ac- 
countinbg for the statistical errors in the am measurements, 
appears to be < 5%, independently of the tracer. 

Such a small scatter in the am — M200 relation poten- 
tially makes am a very good proxy for the mass, via inver- 
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Table 2. Masses (in 10 15 Mq) of clusters at z=0 predicted from 
the a\£i — M200 relation for two values of trio, 300 m/s and 900 
m/s. 



sion of cq. Q: M 20 o/10 15 M Q = (am/ Am)^^ /h(z). How- 
ever, Aid and a are significantly different for the different 
tracers (DM particles, substructures, galaxies). Using the 
values obtained for one tracer to infer cluster masses from 
the (Tid of a different tracer can lead to systematic errors 
of up to ~ 30%. In comparison, the effect of using different 
baryonic physics for the same tracer has a much smaller ef- 
fect on the mass estimates obtained from o\d (< 7%, see 
Table [2]). The presence of scatter, even though small, makes 
possible the applicability of the scaling relation only in a 
statistical sense, not for mass estimates of single objects. 

The results of this paper show that good knowledge of 
the (tid — M200 relation in 6D phase space is fundamental 
before one could apply this relation to observational sam- 
ples. Projection effects and the presence of interlopers can 
significantly a ffect the accuracy and reliabilit y of the mass 
estimate (e.g. ICenl Il997l ; iBiviano et"alll2006h . We plan to 
consider these effects in detail, in a forthcoming work, us- 
ing simulations with a proper implementation of baryonic 
physics and galaxies as tracers. 
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Figure Al. Same as Fig.[3]for the DM set. 
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Figure A3. Same as Fig.[3]for the CSF set. 
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Figure A2. Same as Fig.[3]for the NR set. 

APPENDIX A: PLOTS FOR OTHER MODELS 

In this appendix we report the plots showing the redshift 
dependence of the slope and the intercept of eq. (fTJ) as well 
as the velocity distributions of the tracers for all simulation 
sets not already shown before in the main text. 
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Figure A4. Same as Fig. [9] for the DM set. 
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Figure A5. Same as Fig. |9]for the NR set. 
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